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Abstract 

We have calculated the decay constants of B and D mesons with lattice 
QCD. We use an O(a) improved action that takes light quark actions as 
a starting point, tuned so that it can be directly applied at the physical 
masses of the b and c quarks. Our results are Jb = 164 ~t\f ±8 MeV, /s s = 

185 i 1 ! ± 9 MeV, f D = 194 ±$ ± 10 MeV, and f Ds = 213 ±{f± 11 MeV 
in the quenched approximation. The first error in each case is statistical, 
and the second is from perturbation theory. We show that discretization 
errors are under control in our approach, and smaller than our statistical 
errors. The effects of the quenched approximation may raise our quenched 
result by up to 10%. 



1 Introduction 

This paper uses lattice QCD to calculate the decay constants of B and D mesons 
using an 0(a) improved action that takes light quark actions as a starting point, 
tuned so that it can be directly applied at the physical masses of the b and c 
quarks. The B meson decay constant, fs, is of particular phenomenological 
interest, since it is responsible for most of the uncertainty in current determina- 
tions of the Cabibbo-Kobayashi-Maskawa matrix element Vtd from BB mixing. 
Recent calculations of the decay constants of B mesons give much lower values 
than the earliest calculations done in the static approximation because of sev- 
eral effects, all of which happened to be negative. Perturbative corrections turn 
out to be large and negative. The static approximation has larger statistical 
errors than methods with propagating quarks. This led to contamination from 
excited states, which raised the estimate of fg. The use of smeared sources 
and propagating quark methods mitigate this problem. The 0(a) finite lattice 
spacing corrections likewise turned out to be large and negative. 

Conversely, the earliest results for lighter quarks tended to be too low. They 
were done using naive light quark methods with an incorrect and singular quark 
field normalization which forced Jm^/M - * in the heavy quark limit. 
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The combination of a too high /bVMb from the static approximation and 
a too low fDy/Mo from light quark methods led to very large estimates of 
the 1/M corrections to /mV~M in the static limit. With more recent results, 
including those presented here, the 1/M corrections are much reduced. Refs. [jlj 
and Q contain reviews of some of the early work. 

Currently convenient inverse lattice spacings are no larger than a -1 w 2 or 
3 GeV. Therefore, discretization errors that go like the quark mass in lattice 
units, ma, to a power are unpleasantly large for the c quark and completely out 
of control for the b quark. This means that standard light quark formulations 
for lattice fermion actions cannot be used unaltered for the b quark. 

There are several ways of approximating heavy quarks in lattice QCD cal- 
culations with control over discretization errors. These include Nonrelativistic 
QCD §, |, | (NRQCD), the static approximation [?] , and the approach of 
this paper, described in Refs. |8j |9j. The methods vary significantly both in 
ease of application and in suitability for various calculations. The last takes the 
light quark actions of Wilson jl0| and Sheikholeslami and Wohlert (SW) jll| as 
its leading approximation but adds correction operators that end up resembling 
those of NRQCD rather than those of the standard Symanzik improvement 
program. Sec. || contains a general discussion of the various methods for heavy 
quarks. It comments on our bettered understanding about which methods work 
best in which situations, not all of which was expected in advance. Sec. [| treats 
in detail the particular form of the action used in this paper. Our numerical 
results are presented in Sec. 0, and compared with other recent results in Sec. Ex. 



2 Methods for Heavy Quarks on the Lattice 

This work uses a formalism for propagating heavy quarks that reduces to the 
ordinary light quark formalism in the light quark limit ^ . We begin by com- 
paring the various lattice approaches to heavy quarks. NRQCD is based on an 
expansion in nonrelativistic operators (rotationally invariant but not Lorentz 
invariant) similar to that used in calculating relativistic corrections in the hy- 
drogen atom. It can be thought of as arising from a discretization of the ac- 
tion arising after a Foldy-Wouthuysen-Tani (FWT) transformation of the quark 
fields: 

i/> -> exp(0A7i)V>, (1) 
with 8 chosen so that 



£>o7o + m 
E>2 _ (p2)2 

2m 8m 3 



(2) 



The rest mass term does not affect the dynamics of nonrelativistic particles and 
is conventionally removed. Increasing accuracy is achieved by truncating the 
series with increasing numbers of terms. 

In B physics, one can use the simplest of all the methods, the static approx- 
imation, which is the truncation of the preceding series to a single term. Then 
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the heavy quark propagator is a simple Wilson line in the time direction. It 
is most useful for the heaviest quarks. It is not much used recently because it 
has a much worse signal to noise ratio than methods using propagating quarks, 
which is clear in retrospect but was not foreseen. 

The third method, used in this paper, can be thought of as arising from a 
partial FWT transformation: 

i/> exp(0'A7i)^ (3) 

and 

fl+m — > L»o7o + m + ctiATi 
D 2 



where 8' < 9. This appears to be a crazy thing to do, producing an action that 
combines the defects of the transformed and the untransformed actions. On the 
other hand, it turns out that this is the action we have been using for a long 
time. The Wilson and SW actions have just this form. In addition to the usual 
iplpip term, the i)jD 2 ip term added to cure the doubling problem also contributes 
to the kinetic energy, as in Eqn. ^| The relative strengths of the two terms 
change as ma is increased. 

For the Wilson action, oi = 1, and a 2 — ma. These have the property that 
the nonrelativistic ipD 2 ip term takes over automatically from the Dirac-style 
kinetic energy term vpflip as ma — > oo in the Wilson action. The Wilson action 
turns into a nonrelativistic action in the large mass limit. The heavy quarks in 
heavy-light mesons are highly nonrelativistic. (p/m c h ~ Aqcd/wic/i ~ -2 and 
p/mb ~ .06.) Therefore, it is the Wilson term, rather than the Dirac term, that 
contributes most to the heavy quark kinetic energy. 

We can write a lattice energy-momentum relation 

E 2 = M? + ^p> + ..., (5) 

where the "rest mass" 

Mi - E(0), (6) 
and the "kinetic mass" 

M 2 " x = {d 2 E/dpf) p=0 . (7) 

The Wilson action and most other actions have the property that M\ does 
not equal M 2 for ma ^ 0. For nonrelativistic particles, the rest mass does 
not affect the dynamics and the kinetic mass governs the leading important 
term. Therefore, the rest mass is normally simply omitted from the action in 
NRQCD and in the static approximation, although there is no harm (and no 
benefit) in including it. However, it cannot be set to zero if one wishes to 
recover a sensible light quark limit. It is easy to find a Wilson style action 
which does satisfy M\ = M 2 by letting the hopping parameter for the time 
direction, « t) differ from the one in the spatial directions, k s . The two k's can 
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be separately tuned so that Mi = M2 with no loss of predictive power, by 
requiring Lorentz invariance in Eqn. |5|. Nevertheless, since Mi doesn't affect 
the physics of nonrelativistic systems, the Wilson and SW actions can be used 
without problems for nonrelativistic systems as long as the kinetic mass, and 
not the rest mass, is used to set the quark mass. We can use Kt = ft s with an 
incorrect Mi and a correct Mi- 

The Wilson and SW actions and the action of NRQCD all use one-hop time 
derivatives. When the action includes a large rest mass, this is a requirement, 
since two-hop time derivatives introduce new states with complex energies when 
ma > 1. Therefore, heavy-quark improvements of the Wilson action cannot 
follow the conventional Symanzik program of adding two-hop space-like and 
time-like interactions at 0(a 2 ), but must follow NRQCD and correct only the 
spatial interactions. The existence of the transfer matrix and the Hamiltonian 
ensures that this is possible. 

The parameters of the action in our approach must have nontrivial mass 
dependence, just as those of NRQCD do. For ma > 1, the wave function nor- 
malization, the relation between the physical mass and the hopping parameter, 
etc., are completely different from their ma — values. When ma < 1, this 
mass dependence may be expanded in power series. For the Wilson and SW 
Lagrangians, this yields just the usual series of operators, with the same co- 
efficients. If the mass dependent coefficients in our style of interpreting the 
Lagrangian, 

C = m(am)ipip + z{am)i\)Jp^ + c{am)ipa ^F^ip + . . . , (8) 

are expanded, the usual series of powers of ma multiplying identical operators 
must result: 

C = miptp + aim 2 ''a^pip + . . . 

+ ziplpip + zimaiprfhp + Z2m 2 a 2 ipl/)ip + . . . 

+ cswipcr^F^ip + c 2 ma4><J lllJ Fi lv ip + (9) 

Sometimes in the Symanzik program, i/jCT^F^ip and maipa^F^t/j are spoken 
of as if they were totally unrelated operators. Even for light quarks, it makes 
more sense to think of them as different terms in an ma expansion of an ma 
dependent coefficient analogous to the expansion in a s of an a s dependent 
coefficient. The fact that coefficients in the usual approach blow up as ma — ► 00 
is a property of the expansion and not of the required functions in Eqn. || 
themselves. They stay well behaved if sensible normalization conditions are 
applied. 

With the added ingredient of decoupling the timelike and spacelike parts 
of operators, it becomes possible to formulate an action that is systematically 
improvable, even for large ma. While for ma > 1 the action becomes very similar 
to NRQCD in its behavior, for ma < 1 it may be regarded as a resummation of 
the usual operators of the Wilson and SW actions to all orders in m. 

For physics involving the b quark, NRQCD methods are often easier and 
therefore more accurate. Taking the ma <C 1 limit is not possible, but cor- 
rection operators are not too hard to organize and add. This is particularly 
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important in the T system, where momenta are rather large, and correction 
operators are important in obtaining such things as the correct spin-dependent 
spectrum. In B mesons, on the other hand, the b quarks are extremely nonrel- 
ativistic. 0/c) 2 - (0.3 GeV/5.0 GeV) 2 - 0.3%. Therefore, only the first few 
operators in the action, which are the same in the two methods, are important. 
As long as these are normalized in the same way, it matters little what addi- 
tional operators are hanging around. Therefore, in B mesons the effects of the 
differences between NRQCD and our heavy quark methods are minimal, and 
the two methods should yield nearly identical results. 

For physics involving the c quark, heavy quark methods that can recover a 
relativistic form of the action are often easier and more accurate. The series 
in v 2 is less convergent, so eliminating the need for it by taking a toward zero 
while recovering a relativistic action may be more convenient. For hadrons 
containing charmed quarks, it is possible to do calculations with the Wilson 
and SW actions even with the old interpretation of the coefficients. However, 
since m c ha ~ bah-QCD, our ability to sum up the required series in ma exactly 
is likely to produce a faster approach to the continuum limit than naive light 
quark methods. 

For physics involving cither the b or c quarks, both NRQCD and our method 
can be used successfully. Even when one works better or worse than the other, 
we still learn something. 

3 The Action Used in this Calculation 

We have used the approach for heavy quarks outlined above to calculate the 
decay constants of the D, D Sl B, and B s mesons. We start with an action 
corrected to 0(a), which in general is 

S = ^V'nV'n 
n 

- ~ 7o)£4,oV'„+o + </Wo(l + 7o)f^, ^n 

n 

n;i,j,k 

+ ic E K s ^2'ip n a iE n . ii tp n . (10) 

n;z 

The second and third terms are the timelikc and spacelike pieces of the kinetic 
energy of the Wilson action. The fourth and fifth terms are the spacelikc and 
timelike pieces of the SW correction operator. 

The action used to calculate /mVM may be simplified from this form. As 
previously noted, when ma ^ 0, we find that Mi ^ M 2 and that higher orders 



5 



in p in Eqn. || do not vanish as they should. Explicitly, we define 

ma =-L-[l + 3C], (11) 

lK t 

where ( = k, 8 / Kt- Then we can calculate the free lattice propagator and expand 
in pa to derive 

M x a = log(l + ma), (12) 

and 

1 2C 2 C 



M2Q, ma(2 + ma) 1 + ma 



(13) 



For nonrelativistic fermions, the physics of the system is influenced hardly at all 
by the rest mass Mi. M2 controls the leading order dynamics. We can therefore 
set 

Kt = K s = K, (14) 

but we must tune k so that M2 equals the physical quark mass. 

Likewise, the physics of heavy-light mesons is insensitive to the final, spin- 
orbit term in Eqn. [n]. (This not true of some other quantities in heavy quark 
physics, such as the fine structure of quarkonia.) We can therefore set 

c-e — cb = c, (15) 

and tune c to get the physics of the fourth operator in Eqn. [TJj, the cr-B operator, 
correct. At tree level, this requires the same value that it has for light quarks 
in the SW action, c = csw — 1- 

The decay constant /m parametrizes the matrix element 

(0\A fl \M,p)=p fl f M (16) 

of the axial vector current between the pseudoscalar meson M and the vacuum. 
Here the state \M, p) has the standard, relativistic normalization. Like the 
action, the axial vector current must be specified through 0(a). At tree 
level we take 

= 2^K t j l Kt,iZ A ^hl t ,l5^h (17) 

where the subscripts h and I denote the heavy and light quarks, is a mass- 
dependent (re)normalization factor, and ^ denotes a (rotated) field, specified 
below. 

The factors \/2k arise from naive quark wave function mass dependence in 
long historical use. Wilson used the field normalization 

nai' 
z wf 



= V2k (18) 



for massless fermions [|o| . It is easy to see that this is correct only in the ma = 
limit that Wilson was considering. Away from ma — 0, a straightforward 
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examination of the free propagator shows that, if the quarks are to have their 
conventional canonical normalizations, the correct normalization is 



z wf = VI -6k (19) 
for the Wilson and SW actions, or 



z wl = VI ~Gk s (20) 
if K t ^ K s . These are equal to the conventional but incorrect normalization 



for ma = 0: V2k = vl — 6k = 1/2. As mha — ► oo (nt,h — > 0), the naive 
normalizations force /m incorrectly to zero, unless the mass dependence of 
compensates accordingly. Indeed, the early estimates of fu that were too low 
stem from the naive assumption Za ~ 1, even when m^a ~ 1. Combining 
Eqns. [l^ and we see that 



^ = \ Hi^H 2^77- (21) 



= ^/(l + m h a)(l+mia), (22) 

where the second expression follows from Eqn. 

The remaining ingredient of Eqn. ^ is the rotated field 

* = (l + di 7 ,A)^ (23) 
where i is summed over the three spatial directions and, at tree level, 

qi + ma) 1 (24) 

The second term of di is the same as the 1/m correction present in NRQCD. It 
is a consequence of the Foldy-Wouthuysen-Tani transformation in Eqn. [I]. The 
first term of d\ allows for the fact that the implicit Foldy-Wouthuysen-Tani 
transformation of Wilson-like actions takes 0' 7^ 9. The numerical value of d\ 
depends on the tuning of C, = K s /Kf In most cases, and in particular when 
C = 1, the contribution to Jm proportional to di is 0(mka 2 ) when ma -C 1, 
but 0(A/m) when ma ^S> 1. In the latter region it is essential to include it, if 
one is interested in the 1/mh corrections to the static limit. It also makes sense 
to include the rotation when ma < 1. 

The dominant mass-dependence exhibited in Eqn. ^l] persists beyond tree 
level ||. Quantum effects require further terms to be added to the right-hand 
side of Eqn. [TJ] to construct an C(a)-improved axial vector current. In our nu- 
merical work we neglect them and we set di to its tree-level, mean-field improved 
value. 

Thus, to O(o°) and 0(a x ), we can use precisely the Wilson and the SW 
actions for much heavy quark physics. In this paper we will use the SW action: 

S = y^^ni>n 
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+ ^CsWK^TpnV^K^v^n- (25) 

There are two key differences between our use of this action for heavy quarks 
and the naive approach. 

1. The use of the correctly normalized quark fields and currents, and the 
additional three dimensional field rotation for the heavy quark in the cur- 
rents. 

2. The use of the kinetic mass rather than the rest mass to set the quark 
mass. 

At higher orders in a, it is no longer true that ordinary forms of the light quark 
action may be used for heavy quarks. Two-hop corrections for time derivatives 
cannot be used since they introduce new states with energies which become 
complex when ma > 1. 0(a 2 ) corrections must have a nonrelativistic form as 
in NRQCD. 

Having determined the coefficients in the action at tree level, we now dis- 
cuss their renormalization. In general, the coefficients c„ of operators O n in a 
Wilsonian effective action depend on the physical parameters of the theory and 
on the renormalization scale. For lattice QCD we can write 

S{m/A QC D, Aqcd a ) = X! c ™( ma ' a s)0„, (26) 

n 

where m and a s are the bare parameters of the lattice theory. In the usual 
Symanzik approach, these coefficients are expanded both in a s and in ma. It 
is usually an advantage if an expansion in a s can be avoided. Likewise, as 
remarked earlier, there is an advantage in not expanding in ma. The expansion 
in ma is the source of the breakdown of the method for ma > 1. If the c n 
are not expanded in ma, applying standard normalization conditions yields well 
behaved c n for large ma, as shown at tree level for a few important cases above. 
It is of course necessary that this is true for the c„ in the full quantum theory 
as well. 

Most quantum corrections done so far have only been calculated in per- 
turbation theory. It is clear that the tadpoles, which dominate many lattice 
perturbation theory calculations, create well-behaved contributions to the renor- 
malization of these actions at all quark masses. Tadpole-improved perturbation 



theory 12 suggests that the dominant perturbative correction to the renormal- 
ization of coefficients of operators is given by short distance tadpoles. They can 
be estimated by noting that the link = 1 + igoA^ — ^g^A^ + . . . does not 
fluctuate around 1, as is implicitly assumed in the usual perturbative expansion, 
but around some smaller value, uq- The so-called mean link uq can be estimated 
from the expectation value of any link-containing operator such as the plaquette 
((Up) = Uq) or the expectation value of Un itself in Landau gauge. By count- 
ing links in the actions Eqn. O or Eqn. |23 one can see that the usual tadpole 
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Figure 1: The current normalization factor 2nZy calculated nonperturbatively 
from the quantity ( J/tp\ V4 1 J/ip), where the (unnormalized) local charge V4 = 
ipjoip should simply count the numbers of fermions when properly normalized. 
The upper curve is 1 — 6k, from mean-field improving Eqn. [l^. The lower curve 
is the naive ansatz 2k. 



estimates of corrections to operator coefficients k — kuq and csw = csw u o are 
still valid for heavy quarks even though the one-loop coefficients differ slightly. 

The one-loop coefficients which have been calculated so far have the required 
property that they stay small and well-behaved as ma is varied from to 00 [ p"3| . 

The part of the one-loop correction to the local heavy- light axial vector 
current arising from the leading operators in the SW action has been calculated 
by Aoki et al. (IB). It also is well-behaved for all quark masses, and has been 
incorporated into our results. The part of the correction arising from the 0(a) 
correction to the axial vector current is currently being calculated |L6| . 

One important quantity can be easily computed nonperturbatively, to see 
that the required good behavior persists to all orders of perturbation theory 
This is the quark wave function renormalization, which is required for all current 
renormalizations. This is part of the vector current, which is easy to normalize 
nonperturbatively by inserting the charge into any physical state, such as a IS 
quarkonium state g. Fig. [I] shows the result of such a calculation, compared 
with the tadpole improved tree-level correct normalization and the naive nor- 
malization. The vector current normalization factor Zy is defined analogously 
to the axial vector normalization factor Za in Eqn. It is easy to see that 
the good behavior of the correct normalization in Eqn. [19| is preserved when full 
quantum effects are included nonperturbatively. 
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4 Numerical Results 





G.l 


5.9 


5.7 


Volume 


24 x 48 


16 x 32 


12 x 24 


(!nn n Q~nra t 1 on ^ 


100 


350 


300 


csw 


1.46 


1.50 


1.57 


Qv(l/a) 


0.222 


0.259 


0.330 


a\r(2la\ 


0.171 


0.192 


0.227 


ay(7r/a) 


0.149 


0.164 


0.189 


u 


0.8816 


0.8734 


0.8608 




0.1372 


0.1382 


0.1405 




0.1373 


0.1385 


0.1410 




0.1376 


0.1388 


0.1415 




0.1379 


0.1391 
0.1394 


0.1419 




0.050 


0.050 


0.050 




0.099 


0.093 


0.070 




0.126 


0.1227 
0.126 


0.089 
0.110 
0.119 


a- 1 (/k) (GeV) 


2.21 ±1 


1.571 1 


1.01 1 i 


a~ l {lP- IS) (GcV) 


2.62 ±1 




1.1611 



Table 1: Parameters used in the numerical calculations. 



We turn now to a discussion of our numerical results. A preliminary version 
of our results appeared in Ref. fL8|| , In Table [l] we show some details of our 
calculations. The light quark propagators are the same ones used to determine 
the light quark masses in Ref 0| . They were calculated with the SW action 
using a mean-field improved coefficient csw- cty(l/a), av(2/a), and ay(7r/a), 
defined as in Ref. |T^j , are values of the strong coupling constant at several 
relevant scales, uo is the mean link used in tadpole improvement, ki and 
are the hopping parameters used for the light and heavy quarks, respectively. 
a -1 (fx) and a _1 (lP-lS) are the lattice spacings in physical units as determined 
by fx and the IP-IS splitting of charmonium, respectively. For (3 — 6.1, 5.9, 
and 5.7, the hopping parameters closest to the b hopping parameters are k& « 
0.099, 0.093, and 0.089 respectively. Those closest to charm are n c h « 0.126, 
0.1227, and 0.119 respectively. The hopping parameters closest to strange are 
0.1373, 0.1385, and 0.1405 respectively. Only modest interpolations in quark 
mass are necessary for B s and D s mesons. B and D physics must be obtained by 
chiral extrapolation to the physical light mass limit. Charm and bottom k's were 
obtained by demanding that spin- averaged kinetic masses of the IS charmonium 
and bottomium states match experiment. For the strange quark, the hopping 
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parameter was fixed from the kaon mass. We have also done calculations with 
static propagators for the heavy quarks, which corresponds in our notation to 
k s = and Kt = 1/2. 



1.0 



0.8 



0.6 



cfl 0.4 ~ 



0.2 



(3=5.7 
(3=5.9 
(3=6.1 



o.o 



0.01 



ma: 



0.03 
=ln( 1- 



1/2k 



0.05 
-1/2k 



0.07 



Figure 2: Chiral extrapolations of JbVMb to m 
(center), and (3 — 6.1 (bottom). 



at P = 5.7 (top), = 5.9 



Masses and matrix elements were obtained from minimizing x 2 's using the 
the full correlation matrix. Decay constants were obtained by dividing matrix 
elements by the square root of the experimentally measured meson mass. Sta- 
tistical errors were calculated with the bootstrap method, using 1000 bootstrap 
samples in each fit. Chiral extrapolations of /bV^b at (3 = 6.1, 5.9, and 5.7 
are shown in Fig. ^[ No problems with \ 2 using linear fits were observed with 
the k;'s shown. In runs at lower quark masses, which we did not use here, 
occasional exceptional configurations and poor x 2 's began to appear. 

The possibility of contamination of the results by excited states was checked 
by testing for consistency of one-, two-, and three-state fits using S function, IS, 
and 2S sources at both source and sink. In the sources, one quark is smeared 
with the supposed wave function and the other is a delta function. Shapes of 
the IS and 2S sources were taken from lattice Coulomb gauge valence quark 
wave functions of the B and D mesons. Fig. ^ shows a comparison of one- 
and two-state fits. The measured energy splitting between the ground and first 
excited states and the observed agreement between the ground state energies 
obtained in one- and two-state fits suggests that the contamination from excited 
states is small in the region used for fitting. The measured values of Jm^/M iri 
one-state fits and two-state fits agreed to within statistical errors. 
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Figure 3: Comparison of results of one- and two-state fits at j3 = 5.9, Kh = 
0.1227 (near K c h), and m = 0.1385 (near k s ). Good consistency between the 
two is observed in the fit range for the lightest state. 



The dominant finite lattice spacing errors for the Wilson action are 0(a). 
When these have been removed, as they have been in the SW action, the re- 
maining finite a dependence arises from a combination of 0{a 2 s ), 0(a s a), and 
0(a 2 ). This means that attempts to extrapolate away any remaining finite a 
dependence are uncertain, since we do not have a theory of the functional form 
of the remaining dependence. The ideal situation when using improved actions 
is negligible dependence of the results on a. As shown in Figs. |^ and g, to within 
statistical errors, this has been achieved if fx is used to determine the lattice 
spacing. The results extrapolated to zero lattice spacing are consistent with the 
results at the smallest lattice spacings. Thus, within our statistical errors we 
do not find any evidence for discretization errors. 

The least understood source of uncertainty is due to the quenched approxi- 
mation (the omission of sea quarks) . One effect of the quenched approximation 
is that different answers for decay constants will be obtained depending on which 
physical quantity is used to set the lattice spacing in physical units. The ap- 
pendix tabulates results using fx, fw and the IP-IS splitting in charmonium 
to set the scale. The well-measured physical quantity which most resembles 
the heavy-light decay constants is fx- It therefore offers the best chance that 
some statistical, systematic, and especially quenching errors in fu will cancel 
out. It is relatively easy to determine numerically. (/„., for example, requires 
an additional chiral extrapolation which introduces larger statistical error and 
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Final Value \ 2 /d.o.f 



!b 


164 ± If ±8 MeV 


1.32/1 


f Bs 


185± 1 | =b 9 MeV 


2.06/1 


fn 


194 ±$±10 MeV 


0.77/1 


fn. 


213 ±\j ±11 MeV 


1.01/1 


/d/Ib 


1.19+ § 




Id s 1 Jb 3 


1.14+ ! 




SbJSb 


1.1311 




Jd 3 1 /d 


1.10+ ! 




UbJSbVUdJSd) 


1.031 j 





Table 2: Final results for the decay constants and ratios in the quenched approx- 
imation. The first error for the decay constants is from statistics, the second is 
systematics, mostly from perturbation theory. An additional uncertainty of per- 
haps 10% is present in the decay constants due to the quenched approximation, 
as discussed in the text. The error in decay constant ratios is statistical. 



unreliability. Our K consists of two degenerate quarks of mass (m s + md)/2,) 
Taking a -1 from fx has the pragmatic defect that it yields /m's that lie at the 
bottom of the range given by all of the standard methods of setting the lattice 
spacing, so that removing the quenched approximation is more likely to move 
the results up than down. One way of estimating the effect of the quenched 
approximation is to examine the spread of values obtained by setting the lattice 
spacing in various reasonable ways. However, as shown in the appendix, quanti- 
ties similar to /m tend to give similar results, while quantities more dissimilar, 
such as quarkonia spectra, can give results quite a bit different. It is not easy 
to quantify what is "reasonable" . The most quantitative estimate so far of esti- 
mating the effects of quenching with actual unquenched calculations was given 
in Ref . |l9| . It compares results for quenched and unquenched calculations done 
in similar ways. They report quenching uncertainties of around ±q% for the D 
and D B mesons, and ±J°% for the B and B s . 

We now summarize our results and uncertainties. In Figs. || and |E| we show 
results for /dV®d, Id s \J Md s , Sb\[Mb, and /b„ \J Mb b , using fx to set the 
lattice spacing in GeV. Our final results for the decay constants and their ratios 
are shown in Table ^. Statistical errors in the decay constants are obtained from 
bootstraps over extrapolations to a = linearly in a. The x 2 ' s refer to these 
extrapolations. Statistical errors in the ratios are obtained from bootstraps over 
the ratios of extrapolations thus obtained. 

A common way of estimating the uncertainties from perturbation theory is 
to do the perturbation theory at a range of plausible scales, say 1/a to ir/a. 
That method does not give a sensible result in this case. The scale variation 
has a much larger effect at large a than at small a. When the results are 
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Figure 4: /bVMb and /dV^d extrapolated to zero lattice spacing. The 
statistical error bands are 1 a ranges of linear extrapolations found in bootstrap 
runs. 
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Figure 5: Jb 3 \JMb s and fn s \/Md~ extrapolated to zero lattice spacing. 
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extrapolated linearly to a = 0, they give the same answers to within a per cent, 
even though the perturbative effect being extrapolated is not at all linear in 
a. This is clearly an underestimate of the perturbative uncertainties. These 
are due to 0(a^) terms of order 3% plus the remaining 0(a s ) term which is 
currently being calculated. Assuming that this term has roughly the same value 
that it has in NRQCD, we have taken 5% for the perturbative uncertainty. 

The remaining uncertainties in the quenched approximation are likely to 
be smaller than the statistical and perturbative uncertainties. We do not see 
evidence in Figs. f| and || for finite lattice spacing effects within our statistical 
errors. We estimate that the effects of uncertainties in the determination of the 
heavy quark mass are about a per cent or two. Finite volume errors are likewise 
expected to be small for B and D mesons on our lattice volumes. 

Quenching uncertainty estimates of around 10% are common, based on 
ranges of results using different quantities to set the lattice spacing. The 
crude, but quantitative estimates of quenching effects in Ref. fl9|| are within 
the range of the less quantitative estimates. If we assume that the effects of 
the quenched approximation are similar for different quark actions, and apply 
the largest of their quenching uncertainties to our quenched results, we obtain, 
f B = I64±\j±8± 16 MeV, f Bt = 185^1 ±9i^ MeV, f D = 194±i3±10t 2 g 

MeV, and f Bs — 213+if ± lli 2 o MeV, where the uncertainties are statistical, 
systematic other than quenching, and quenching. The last uncertainty is less 
quantitative than the first two. 

Perturbative errors cancel in ratios of decay constants. It is possible, but 
not proven, that deviations of ratios from one have quenching uncertainties of 
the usual 10-20% size. In that case they would be significantly smaller than 
our statistical uncertainties in these quantities. We have therefore quoted only 
statistical errors for these quantities. 

Our value for f Ba is compatible with the world average experimental value 
f Ds = 241 ±21 ±30 MeV [§0|. The double ratio (f B J f B )/(f D J Id) is within a 
few per cent of one, in accordance with expectations from chiral symmetry and 
heavy quark symmetry pl[ . 

5 Comparison with Other Results 

Refs. [^2[ ^||, |24|, [25) contain recent reviews of heavy-light decay constants. Most 
recent calculations are consistent, with values of around 160 MeV for f B and 
around 195 MeV for f B . The JLQCD collaboration has completed a calculation 
of heavy-light decay constants using methods almost identical to ours (26) . They 
have used the p mass rather than fx to set the lattice spacing, but that makes 
a negligible difference. The results should agree to within statistical errors, 
and they do. They obtain f B — 163 ± 9 ± 8 ± 11 MeV, where the errors arise 
from statistics, systematics other than quenching, and quenching. These results 
are also reasonably compatible with results using Wilson fermions without the 
0(a) correction which rely on extrapolating linearly to a = to remove the 
leading errors. The MILC collaboration has obtained in preliminary results 
f B = 153 ± 10 i 13 i X q MeV, where the errors arise from statistics, systematics 
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other than quenching, and quenching fL9| , Allton et al. obtain higher results 
than the others, with f B = 180 ± 32 and f D = 221 ± 17 MeV, though not 
inconsistent within errors |27 . They have used a tree level coefficient for the 



clover term which is 2/3 the value after quantum corrections. This should leave 
some residual discretization error, but they do not try to extrapolate it away, 
and it is responsible for their large error estimate for fs- In addition, they 
have not performed the fs calculation at the physical b quark mass, but have 
extrapolated from the region m\,a < 1. This may add a source of uncertainty 
which is hard to quantify. 
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Figure 6: Comparison of the continuum limit results of this paper (filled sym- 
bols) with the results of NRQCD (open symbols). The errors shown are statis- 
tical. 



As remarked earlier, for the B meson system, the methods of NRQCD are 
also very similar to ours. This is less true for the D, ip, and T systems. It is 
therefore interesting to compare recent results from the two methods. Fig. [| 
shows results for JmV^M as a function of 1/M from this paper and from a recent 
NRQCD calculation (28). The NRQCD calculation includes a full O(a) correc- 



tion and a full one- loop correction. This paper uses a full O(a) correction, but 
is missing the part of the one- loop correction arising from the 0(a) correction to 
the current, which has not yet been calculated. The one- loop corrected results 
are close to each other in the region of the B mass, where we have reported 
our main results. However, the NRQCD one-loop corrections increase greatly 



1G 



for larger M, becoming as large as —25% in the static limit, the result of a 
surprisingly large contribution from the 0(a) current correction. 

For all of the heavy-light decay constants, the largest source of uncertainty 
uncertainty, and the one under least quantitative control, is the quenched ap- 
proximation. The uncertainties from all other sources of error are now less than 
the approximately 10% range usually associated with the quenched approxima- 
tion. 
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A Appendix 

We present in Tables |^ - [n] results for the decay constants obtained with various 
fitting methods. We give two tables for each of the four decay constants. The 
first table contains results for each of the three lattice spacings. The second table 
gives the results of extrapolation to the continuum limit, along with the \ 2 of the 
fit. In each case, we performed the analysis with quark masses determined from 
the spin-averaged IS -onium state, using one of five methods of determining the 
lattice spacing: the IP-IS charmonium splitting, fx from either 2- or 3-state 
fits, or f v from either 2- or 3-state fits. 



(9 = 6.1 = 5.9 = 5.7 



IP-IS 

f n (2-state) 
f„ (3-state) 
fx (2-state) 
fx (3-state) 



o.5ii ±§§ 

0.419 ±H 
0.4171H 
0.394 if | 
0.391 ±?B 



0.455 ii 
0.407i?§ 
0.363 i 2 ^ 
0.371 iH 



0.477ii 
0.440 ±H 
0.384 ±\\ 
0.385 i^ 



0.357i}g 0.387 i 1 ? 



Table 3: /bVMb in GeV'? at = 6.1, 5.9, and 5.7, fitted with the lattice 
spacing set in five different ways. 
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a- 1 f B VM(a = 0) f B (a = 0) X 2 /d.o.f 





IP-IS 


0.490 iff 


213 


2.04/1 


u 


(2-state) 


0.386 ±p 


168 if | 


0.51/1 


u 


(3-state) 


0.349±g| 


152 ±U 


1.63/1 


Ik 


(2-state) 


0.376 ii 


164 ill 


1.32/1 


Ik 


(3-state) 


0.331 ±|| 


144 +« 


3.51/1 



Table 4: /b\/Mb in GeVa and f B in MeV extrapolated to a = 0, fitted with 
the lattice spacing set in five different ways. 



a' 1 (3 = 6.1 (3 = 5.9 (3 = 5.7 

IP-IS 

fn (2-state) 
f v (3-state) 
Jk (2-state) 
Jk (3-state) 

Table 5: f Bs yjM Bs in GeVi at (3 = 6.1, 5.9, and 5.7, fitted with the lattice 
spacing set in five different ways. 



0.582 tf 2 
0.478 t\t 
0.476 ±2? 
0.449 if g 
0.445 ± 3 f 



0.520±^ 
0.465 ±?f 
0.415± 2 | 
0.424± 1 | 
0.408 1 1 



0.545 ii 
0.503 if £ 
0.439 if 1 
0.440 t n 6 
0.442 i^ 



a- 1 f Bs VM(a = 0) / B .(a = 0) xV^-/ 





IP-IS 


0-567 til 


244 ii 


2.80/1 


/it 


(2-state) 


0.444 i^ 


192 tf 


0.59/1 




(3-state) 


0.395 ttl 


170 ifl 


1.81/1 


/if 


(2-state) 


0.430 if! 


185 i 1 ! 


2.06/1 


/* 


(3-state) 


0.374 ±3i 


161 i 1 ! 


6.51/1 



Table 6: f Bs yj M Bs in GcV'a and f Bs in MeV extrapolated to a = 0, fitted 
with the lattice spacing set in five different ways. 
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a- 1 /3 = 6.1 /3 = 5.9 /3 = 5.7 





IP-IS 


0.332 


+ 16 

-23 


0.329±i| 


0.323±i| 




(2-state) 


0.273 


+ 23 

- 9 


0.294l 2 S] 


0.298 i 1 ! 




(3-state) 


0.271 


+ 35 
-24 


0.262± 1 1 


0.261 1 1 ! 




(2-state) 


0.256 


+ 15 

- 8 


0.268 1 f 


0.261 1 1 


/* 


(3-state) 


0.254 


+ 21 
- 7 


0.2581 1 


0.262 1 1 



Table 7: !d\[Md in GeV^ at = 6.1, 5.9, and 5.7, fitted with the lattice 
spacing set in five different ways. 



a- 1 f D VM(a = 0) f D (a = 0) X 2 /d.o.f 





IP-IS 


0.339 ±H 


248 


+ 20 
-24 


0.03/1 


u 


(2-state) 


0.262 1?§ 


191 


+ 28 
-13 


0.43/1 


u 


(3-state) 


0.268±lg 


196 


+ 29 
-15 


0.06/1 


Ik 


(2-state) 


0.265 ±\l 


194 


+ 14 
-10 


0.77/1 


Sk 


(3-state) 


0.249111 


182 


+ 11 
-10 


0.01/1 



Table 8: /d\/Md in GeV^ and fu in MeV extrapolated to a = 0, fitted with 
the lattice spacing set in five different ways. 



= 6.1 = 5.9 = 5.7 



IP-IS 0.380 1^ 0.376111 0.372±i| 

f n (2-state) 0.312± 2 | 0.336 l 2 f 0.343 1 1 ! 

U (3-state) 0.310 ±\% 0.300 1 1 ? 0.300 1 1 ^ 

f K (2-state) 0.293 1 1 ? 0.3061 § 0.301 ± \ 

f K (3-state) 0.290 1 2 | 0.295 ± I 0.3021 | 



Table 9: f Ds v / M Ds in GeV^ at = 6.1, 5.9, and 5.7, fitted with the lattice 
spacing set in five different ways. 
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f Da VM(a = 0) f Ds (a = 0) X 2 /d.o.f 





IP-IS 


0.384 ±H 


274 


+ 22 
-25 


0.01/1 


u 


(2-state) 


0.295 ±H 


210 


+ 28 
-11 


0.42/1 


u 


(3-state) 


0.303±tB 


216 


+ 33 
-14 


0.08/1 


Sk 


(2-state) 


0.299 ill 


213 


+ 14 
-11 


1.01/1 


Sk 


(3-state) 


0.2811H 


201 


+ 12 
- 8 


0.10/1 



Table 10: Jd s \J Mb, in GeVs and fn B in MeV extrapolated to a = 0, fitted 
with the lattice spacing set in five different ways. 
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